ENTROPY SOLUTION THEORY FOR FRACTIONAL 
DEGENERATE CONVECTION-DIFFUSION EQUATIONS 

SIMONE CIFANI AND ESPEN R. JAKOBSEN 

Abstract. We study a class of degenerate convection diffusion equations 
witfi a fractional non-linear diffusion term. This class is a new, but natural, 
_ generalization of local degenerate convection diffusion equations, and include 

^^J , anomalous diffusion equations, fractional conservations laws, fractional Porous 

Cn ■ medium equations, and new fractional degenerate equations as special cases. 

►^ ' We define weak entropy solutions and prove well-posedness under weak reg- 

^ ' ularity assumptions on the solutions, e.g. uniqueness is obtained in the class 

)—^ , of bounded integrable solutions. Then we introduce a new monotone conser- 

^H ■ vative numerical scheme and prove convergence toward the entropy solution 

w~v ' in the class of bounded integrable BV functions. The well-posedness results 

are then extended to non-local terms based on general Levy operators, con- 
nections to some fully non-linear HJB equations are established, and finally, 
some numerical experiments are included to give the reader an idea about the 
P^ ' qualitative behavior of solutions of these new equations. 

< 

■3 

1. Introduction 

In this paper we study well-posedness and approximation of a Cauchy problem 
^^ , for the possibly degenerate non-linear non-local integral partial differential equation 

C^; .^^. idtu + V-f{u)^-{-A)^/^A{u) inQT=M'^x(0,T), 

m ! ^ ' ' |m(x,0) =uo{x) inM'', 

G^ ■ *- 

\l^ I where / = (/i, . . . , /rf) : M ^- R and A : M — > R are Lipschitz continuous with 

lO ' Lipschitz constants Lf and La, A{-) non-decreasing with A{0) — 0, and the non- 

^^ , local operator — (— A)'*'/^ (or g[-] in shorthand notation) is the fractional Laplacian 

defined as 



o 



X 



cd 



-(-A)V2<^(,) = c. P.V. f '^(" + ^-2"/^"''^ dz 

J\z\>0 



|2|d+A 



Vh ' for some constants ca > 0, A G (0,2), and a sufficiently regular function (p. Note 



that A(-) can be strongly degenerate, i.e. it may vanish on a set of positive measure. 
Equation (jl.ip is a fractional degenerate convection diffusion equation, and this 
class of equations has received considerable interest recently thanks to the wide 
variety of applications. They encompass various linear anomalous diffusion equa- 
tions (/ = and A{u) = u), scalar conservation laws [TBI ISSl [311 IMl ISH] {A = 0), 
fractional (or fractal) conservation laws dHH] {A{u) = u), and some (but not all!) 
fractional Porous medium equations [T7] (/ = and A{u) = \u\u^, m > 1), but 
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see also [SI [7] . Equation (|1.1|) is an extension to the fractional diffusion setting of 
the degenerate convection-diffusion equation [51 |SI] 

(1.2) dtu + V ■ fiu)::^AA{u). 

When A{-) is strongly degenerate, equation (|l.ip has never been analyzed before 
as far as we know. 

The literature concerning the type of equations mentioned above is immense. 
We will only give a partial and incomplete survey of some parts we feel are more 
relevant for this paper. For a more complete discussion and many more references, 
we refer the reader to the nice papers pi and [32] . But before we continue, we would 
like to mention actual and potential applications. A large variety of phenomena in 
physics and finance are modeled by linear anomalous diffusion equations, see e.g. 
[m SJ [TT . Fractional conservation laws are generalizations of convection-diffusion 
equations ( (jl.2p with A{u) = u), and appear in some physical models for over-driven 
detonation in gases [12] and semiconductor growth [JT] , and in areas like dislocation 
dynamics, hydrodynamics, and molecular biology, cf. [Tl[3j[T9]. Similar equations, 
but with slightly different non local term, also appear in radiation hydrodynamics 
[57] . Equations like (|1.2p are used to model a vast variety of phenomena, includ- 
ing porous media flow |39] . reservoir simulation [22j . sedimentation processes [B], 
and traffic flow [ID]. Finally, we mention [53] where degenerate elliptic-parabolic 
equations with fractional time derivatives are considered. 

In the non-linear and degenerate setting of (jl.ip . we can not expect to have 
classical solutions and it is well-known that weak solutions are not unique in general. 
In the setting of fractional conservation laws this is proved in e.g. [21 131 133] • To get 
uniqueness we impose extra conditions, called entropy conditions. In this paper we 
will introduce a Kruzkov type entropy formulation for equation p.ip . This type of 
formulation was introduced by Kruzkov in [34j , and used along with a doubling of 
variables device, to obtain general uniqueness results for scalar conservation laws. 
Much later, Carrillo in [5] extended these results to cover second order equations 
like (|1.2p , see also [31] for more general results and a presentation and proof which 
is more like our own. More recently, Alibaud [T] extended the Kruzkov formulation 
and uniqueness result to the fractional setting. He obtained general results for 
fractional conservation laws. In a new work by Karlsen and Ulusoy [32], a unified 
formulation is given that essentially includes the results of Alibaud and Carrillo as 
special cases. In [2 [32] the fractional diffusion is always linear and non-degenerate. 

The entropy formulation we use is an extension of the formulation of Alibaud, 
and it allows us to prove a general L^-contraction and uniqueness result for bounded 
integrable solutions of the initial value problem (|l.ip . Our uniqueness proof relies 
on some new observations and estimates along with ideas from [5] |3I]- From a 
technical point of view, our proof for A G (0, 2) is more related to the conservation 
law (or fractional conservation law) proof than the more technical proof of Carrillo 
for A = 2 (equation p.2p V E.g. we do not need a "weak chain rule" and hence do 
not need to assume any extra a priori regularity on the term A{u). 

In practice to solve (jl.ip we must resort to numerical computations. But since 
the equation is non-linear and degenerate, many numerical methods will fail to 
converge or converge to false (non-entropy) solutions. The solution is to construct 
"good" numerical methods that insure convergence to entropy solutions. In the 
conservation law community, it is well known that monotone, conservative, and 
consistent methods will do the job for you. There is a vast literature on such 
methods, we refer the reader e.g. to |26j and references therein. For non-linear 
fractional equations there exist very few methods and results so far. Dedner and 
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Rhode [IH] introduced a convergent finite volume method for a non-local conserva- 
tion laws from radiation hydrodynamics. Droniou [T^ was the first to define and 
prove convergence for approximations of fractional conservations laws. Karlsen and 
the authors then introduced and proved convergence for Discontinuous Galerkin 
methods for fractional conservation laws and fractional convection-diffusion equa- 
tions in [inilTT]. After that, the authors introduced a convergent spectral vanishing 
viscosity method for fractional conservations laws in [9J. Kuznetzov type error es- 
timates were also obtained in jHlIin]- In this paper, we discretize for the first time 
(jl.l|) in its general form. We introduce a new difference quadrature approxima- 
tion that we prove converges to the entropy solution. The convergence holds for 
bounded integrable BV solutions, and hence we also have existence of solutions in 
this class. Finally, existence of solutions in the wider class of bounded integrable 
function is obtained through approximation via bounded integrable BV solutions 
(cf. Theorem |4J|. 

In many applications, especially in finance, the non-local term is not a fractional 
Laplacian, but rather a Levy type operator g^: 



9i^[(fi\{x) = / (j){x + z)~ (j){x) - z ■ ^ (j){x)l\^\^i d/x(z), 

"'|2|>0 

where the Levy measure /x is a positive Radon measure satisfying 

|zp A 1 ^(dz) < oo. 



'lz|>0 

These operators are the infinitesimal generators of pure jump Levy processes. We 
refer to [H [13] for the theory and applications of such processes and to (55] for a very 
relevant and nice discussion and many more references. The entropy solution theory 
related to such operators is very similar to the one for fractional Laplacians, and 
the first well-posedness results were obtained in [32 . In this paper we extend the 
entropy theory for (|l.ll) to this Levy setting (cf. equation (15. ip ). Our formulation 
is an extension of Alibaud's formulation and is different from the one given in j32j . 
We also treat completely general Levy measures, i.e. our Levy operators are slightly 
more general than the ones in iS? . 

We also discuss the fact that (|l.ip is related to fully non-linear HJB equations, 
see Section [51 We first show an easy extension of results from [3S] : In one space 
dimension the gradient of a viscosity solution of a fractional HJB equation is an en- 
tropy solution of a fractional conservation law. Then we show a new correspondence 
for any space dimension: If u is a viscosity solution of 

then V = g^[u] is the entropy solution of 

vt-gMi'")] =0- 
The relevance of these results are discussed in Section [51 The final part of the paper 
is devoted to numerical simulations to give the reader an idea about the qualitative 
behavior of the solutions of these new equations. 

Here is the content of the paper section by section. The entropy formulation is in- 
troduced and discussed in Section[21 In Section[31 we state and prove L^-contraction 
and uniqueness for entropy solutions of (jl.ll) . The monotone conservative numeri- 
cal method is then introduced and analyzed in Sectional In Section [SJ we extend 
the well-posedness results proved for solutions of (II. ip to a wider class of equations 
where the fractional Laplacian has been replaced by a general Levy operator. In 
Section [S] we show how solutions of equations of the type (jl.ip are related to solu- 
tions of fully non-linear HJB equations, and in the last section, we provide several 
numerical simulations of problems of the form (II. ip . 
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2. Entropy formulation 



In this section we introduce an entropy formulation for the initial value problem 
(ll.ip which generalizes Alibaud's formulation in [J. To this end, let us split the 
non-local operator g into two terms: for each r > 0, we write g[ip\ — gr['~p\ + g''['^] 
where 

r -, / N „ f (P(X + Z) — ip(x) , 

gM{x)^c,P.V. ^^^^Lidz, 

J\z\<r Fl 

rr 1/ N f ip{x + Z) - Lp(x) 

g M(x) = ex r--^^ Az. 

J \z\>r \^\ 

The Cauchy principal value is defined as 

P.V. I '4^{z) Az — lim / v{z) dz. 
Note that, by symmetry. 



P-y- / r-F77TY dz = 



and hence 



^p{x + z) — (p{x) — z ■ \'ip{x) 



(2.1) gA^]{^)^cxP.V. ^ '-^^, ^^dz. 

Whenever (p is smooth enough, the principal value in (|2.ip is well defined by the 
dominated convergence theorem since 

{cx\\Df\\L'^(Bix^r)) /|2|<,. ^JItx dz when A € (0, 1) 1 
f |pVlU~(B(=.,r))4l<,^dz when A e [1,2) J 
The above integrals arc finite because in polar coordinates they are proportional to 

/ -dTX ^'^'^ ^* f°^' ^ ^ (0, 1) and / -— ^ s'^-i ds for A G [1, 2). 
Jo ^ Jo ■^'' 

This estimate also shows that the integral in (|2.ip exists and this leads to an alter- 
native definition of the operator gr avoiding the principal value (i.e. (|2.ip without 
P.V.). This second definition is used e.g. in X. 

Let us introduce the functions r/k{u) — \u — k\, rj'j^lu) — sgn(w — fc), and qk{u) = 
ri'i.{u){f{u) — f{k)) where the sign function is defined as 

1 for s > 

sgn(s) = ■^ for s = 

-1 for s < 0. 

The entropy formulation we use is the following: 

Definition 2.1. A function u is an entropy solution of the initial value problem 
(jl.ip provided that 

i) uei°°(gT)nC([0,T];Li(R'^)); 
ii) for all fc S M, all r > 0, and all nonnegative test functions Lp £ C^{Qt), 



r]k{u)dt(p + qk{u) ■ \7ip + ijAik)^^)) oAf] + ii^iu) g''[A{u)]ip dxdt > 0; 

Hi) u{-, 0) = uq{-) a.e. 

Remark 2.1. By C{[0,T]; L^{M.'^)) we mean the Banach space where the norm is 
given by ||0||c([o,T];Li(R<i)) = maxtg[o,T]{ /nd \H^,t)\ dx}. 
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Remark 2.2. In view of ij and the properties of / and A, 77^(1*), ^^(m), ?7^(j.')(A(u)) € 
L°^{Qt) while A{u) e L°°{QT)nC{[0,T]; L^(R'^)). It immediately follows that the 
local terms in ii) are well-defined. Since 5r[</>] & C!^{Qt) for £ C^{Qt), also 
the (7r-term in ii) is well-defined. Finally we note that g'^ [ip] (x) is well-defined and 
belongs to L°°(R'') for V e £°°(K'*), and to L'^{R'^) for ^ e Li(R'^) by Fubini 
(integrating first w.r.t. x). It follows that g''[A{u)] G L^^iQr) n C([0,r]; ^^(R'')), 
and hence that the g'"-term in ii) is well-defined. 

Remark 2.3. Since u G C{[0,T];L^{R'^)) by i), part iiij implies that the initial 
condition is imposed in the strong L^-sense: 

lim||u(-,i)-wo||Li(Rrf) =0. 

A more traditional approach where initial values u(-, 0) are included in the entropy 
inequality ii) would also work, cf. e.g. [26, Chapter 2]. 

Let us point out that, in the case A € (0, 1) and whenever the entropy solutions 
are sought in the BT^-class, Definition 12.11 can be simplified to the following one: 

Definition 2.2. A function u is an entropy solution of the initial value problem 
(II. ip provided that 

i) ue L°°{Qt) n L°°(0, T; BV{W^)) n C([0, T];L^{W^)); 
ii) for all k €M. and all nonnegative test functions ip G C^{Qt), 

r]k{u)dtip + Qkiu) ■ V(p + ri'f.{u) g[A{u)] (p dxdt > 0; 

Qt 

Hi) u{-, 0) = uo(-) fl-e. 

Note that the non-local term (7[yl(u)] in the integral in ii) is well defined as shown 
in the following lemma. 

Lemma 2.4. If X (z (0, 1), then there is a constant C > such that 



\9[M'^)]\\l^W) < cxCLa\\u 



lLi(Rd)l"'ISV 



Proof. We split the integral in two parts, use Fubini and the estimate 

\u{x + z) - u{x)\ dx < Vd|z||u|By(Rd) 



(cf. Lemma |A.ip . and change to polar coordinates (z = ry for r > and \y\ — 1) 
to find that: 

dxdz < LAVd\u\Bv(R-i) / TTyXT ^^ 



\z\<eJR'' 1^1 J\z\<e Fl 



= LAVd\u\BY^Ts.d) / dSy — = LAVd\u\BV(R'')£^ ^ / ^Sy / 

J|i,|=i Jo ''" J\y\=i Jo 

and 

f f \A{u{x + z)) - Ajuix))] 2La.. ii f ,„ r 

J z>eJR'' 1^1 f^ Jk =1 Jl 



To conclude, we choose e = |ju||Li(R<')|w|j3y/Rdy D 

The following result shows how the two definitions of entropy solutions are in- 
terrelated and how they relate to weak and classical solutions of (II. ip . 

Theorem 2.5. 

i) Definition \2.1\ and Definition \2.'A are equivalent whenever A G (0,1) and 
u G L°°(Qt) n L°°(0, T; BViR"^)) f] C([0, T]; L\R'')). 
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a) Any entropy solution u of (jl.ip is a weak solution: for all ip G C^{Qt), 
udtif + f{u) ■ \/ip + A{u) g[ip\ dxAt = 0. 

Hi) If AG C2(M), then any classical solution u € L°°{Qt) D C{[0,T]; L^R"^)) 
of (ll.ip is an entropy solution. 

Remark 2.6. In Hi) we need additional regularity of A to give a pointwise sense to 
the equation and hence also to define classical solutions. When X € [1, 2) it suffices 
to assume that A G C^, and when A e (0, 1) A G C^ is enough. 

Proof. 

i) Repeated use of the dominated convergence theorem and Lemma 12.41 first 
shows that, when r — > 0, 

gr[Lp\-^Q and g" [Aiu)] -^ g[A{u)] a.e., 

and then combined with this convergence result and Holder's inequality, that Def- 
inition 12.11 implies Definition 12.21 when u is BV. To go the other way, let us note 
that since A{-) is non-decreasing, 

(2.2) sgn(u - k){A{u) - A{k)) = \A{u) ~ A{k)\. 

Thus, if we write 

g[A{u)]=gAA{u)] 

f {A{u{x + z, t)) - A{k)) - {A{u{x, t)) - A{k)) 
+ ^W g^TA d^ 

Je<\z\<r 1^1 

multiply each side by rj'^.{u)Lp and integrate over Qt, we end up with 



V'M 9[A{u)] V dxdi < // <^ r^'kiu) g,[A{u)] ^ 
Jjqt \ 

\A{u{x + z, t)) - A{k)\ - \A{u{x, t)) - A{k)\ 



+ Cxip 

Je<\z\<r 

+ Vkiu)g'^[A{u)]^\ dxdt. 



We now use the change of variables (z, a;) — > (— z, x + z) to pass the test function 
Lp inside the integral e < \z\ < r, and obtain 

^(., t) f \Aiui. + ^,t))-Am-\Aiuix,t))-Am ^^^^^^ 

= // \Aiuix,t))-Aik)\r ^(:i±i4,l^(M dzdxd. 

The entropy inequality in Definition 12. II is finally recovered in the limit as e — > 0. 
ii) Using (|2.2[) and the change of variables {z,x) -> {~z,x + z), 

r/'f,{u{x, t)) g^[A{u{x, t))] (p{x, t) dxdt 



<c,ll pix,t)[ \Mui^ + ^^t))^Mk)\-\Aiuix,t))-Aik)\^^^^^^ 

J\z\>r kl 

\A{u{x,t)) -~Aik)\g''[ipix,t)]dxdt. 

Qt 



DEGENERATE FRACTIONAL ORDER CONVECTION-DIFFUSION EQUATIONS 

Thus, since g = gr + g^ , we have produced the inequahty 

rikiu)dt(p + qkiu) ■ V(p + r]A(k)iA{u)) g[ip] dxdt > 0. 
By this inequahty and the definitions of i] and q, if ±fc > ||ii||/^oo/R), then 



T // {u- k)dt<j> + if{u) - f{k)) ■ V0 + {A{u) - A{k))g[cj,] Axdt > 0. 

By the Divergence theorem and a computation hke in (|2.3p , aU the fc-terms are zero 
and hence m is a weak sohition as defined in ii). 

Hi) Since u solves equation (ll.ip point-wise, for each (x, t) e Qt and aU /c G M, 
we can write 

dtiu-k) + V-ifiu)-f{k))=9,[A{u)] 

f {A{u{x + z, t)) - A{k)) - {A{u{x, tj) ^ A{k)) ^ 

+ CA / r-FTTTT dz 



/ , , lzK+^ 

Je<\z\<r 1^1 

+ g'-[Aiu)]. 

If we muhiply both sides of this equation by ry^ (u) and use (|2.2[) , we obtain 

77^^) at(« - fc) + r7^(«) V • (/(u) - /(fc)) < ^',iu) g,[Aiu)] 

f lAjujx + z,t)) - A(fc)| - \A{u{x,t)) ^ A{k)\ ^ 
+ cW ^p^3^ dz 

Let us now multiply both sides of this inequality by a nonnegative test function ip, 
and integrate over Qt to obtain 

r]kiu) dtip + qk{u) ■ Vifi dxdt 



< II {v'k{u)ge[A{u{x,t))]^ 

\A{u{x + z,t)) - A{k)\ - \A{u{x,t)) - A{k)\ 

|7jdTA 



-Cx(p 

I €<\z\<r 



+ 7j',,{u)g^[A{u{x,t))]^\ dxdt. 

Thanks to (|2.3p . we can pass the test function cp inside the integral e < \z\ < r, and 
so recover the entropy inequality in Definition 12.11 in the limit as e — > 0. D 

3. L^-CONTRACTION AND UNIQUENESS 

We now establish L^-contraction and uniqueness for entropy solutions of the ini- 
tial value problem (|l.ip using the Kruzkov's doubling of variables device [HJ. This 
technique has already been extended to fractional conservation laws (i.e., A{u) = u) 
by Alibaud [1 . The first part of our proof builds on the ideas developed by Alibaud 
(and Kruzkov!), but in the rest of the proof different ideas have to be used in our 
non-linear and possibly degenerate setting. 

Theorem 3.1. Let u and v be two entropy solutions of the initial value problem 
p.l[) with initial data uq and vq. Then, for all t £ {0,T), 

\\u{-,t) -!'(•, i)||Ll(R'') < \\uq -VqW L^j^dy 

Uniqueness for entropy solutions of (jl.ll) immediately follows from the above 
L^-contraction: if uq = wq, then u = v a.e. on Qt- 
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Corollary 3.2. (Uniqueness) There is at most one entropy solution of (11.11) . 

Proof of Theorem \3.1\ 

1) We take u = u{x,t) and v = v{y,s), let xp = ip{x,y,t,s) be a nonnegative 
test function, and denote by r]{u,k), q{u,k), rj'{u,k) the quantities rjk{u), qk{u), 
rj'f,(u). After integrating the entropy inequality for u — u{x, t) with k = v{y, s) over 
(y, s) € Qt, we find that 

T]{u{x, t), v{y, s)) dt^l^ix, y, t, s) 

(3.1) + '?('"(a;, t), v{y, s)) ■ Wa:ip{x, y, t, s) 

+ 7]iA{u{x, t)), A{viy, s))) grl^i; y, t, s)]{x) 

+ T]'{u{x,t),v{y,s)) g''[A{u{-,t))]{x) ip{x,y,t,s) da;didyds > 0. 

Similarly, since rj{u,k) — r]{k,u), q{u,k) = q{k,u), and rj'iu^k) — ~r/'{k,u), inte- 
grating the entropy inequality for v — v{y, s) with k — u{x, t) leads to 

r]{u{x, t), v{y, s)) dsil^ix, y, t, s) 

(3.2) + ^("(a;, t),v{y, s)) ■ Wyipix, y, t, s) 

+ tjiAiuix, t)),A{v{y, s))) grii^ix, ■,t, s)]{y) 

- r]'{u{x, t),v{y, s)) g''[A{v{-, s))]{y) iIj{x, y, t, s) dydsdxdt > 0. 

Let us now introduce the operator 

~rr t \^t \ f ip{x + z,y + z)~ip{x,y) 

J\z\>r 1^1 

Since all the terms in p.ip - (|3.2[) are integrable, we are are free to change the order 
of integration, and hence add up inequalities (j3.ip - p.2p to find that (from now on 
dw = dx dt dy ds) 

(3.3) 

l{u{x, t),v{y, s)) (dt + ds)'4){x, y, t, s) 



+ q{u{x, t),v{y, s)) ■ {W^. + Wy)ilj{x, y, t, s) 

+ r]{A{u{x, t)),A{v{y, s))) gr[ip{-,y, t, s)]{x) 

+ 77(A(u(a;, t)),A{v{y, s))) gr[tlj{x, ■,t, s)]{y) 

+ r^'{u{x,t),v{y,s)) f[A{u{-,t))-A{v{-,s))]{x,y) ij;{x,y,t,s) dw>0. 

In the following we will manipulate the operator g^, while the operators gr will 
simply be carried along to finally vanish in the limit as r — ?► 0. 
Let us use (12. 2p to obtain the (Kato type of) inequality 

r]'iu{x, t),v{y, s)) [(^A{u{x + z, t)) ~ A{v{y + z, s))) - (a(m(x, t)) - A{viy, s))) 
< \A{u{x + z,t)) - A{v{y + z,s))\ - \A{u{x,t)) - A{v{y,s))\, 

which implies that 

(3.4) 

V'iuix, t), v{y, s)) nA{<;t)) - A{v{; s))]{x, y) < f \\A{u{-, t)) - A{v{-, s))^ [x, y). 



DEGENERATE FRACTIONAL ORDER CONVECTION-DIFFUSION EQUATIONS 



Furthermore, we use Fubini's Theorem and the change of variables (z, x, y) 

(— z, X + z,y + z) to see that 



(3.5) 



QtJJQ 



ij{x,y,t,s)f \A{u{-,t))-A{vi-,s))\ ix,y) dw 



\A{u{x,t))-A{u{y,s))\ g''[ij{-,-,t,s)]{x,y) dw. 



iQt JJQt 

To sum up, when used in (I3.3p . (|3.4p - p.5p produce the inequahty 

^{u{x, t), v{y, s)) {dt + ds)ipix, y, t, s) 



(3.6) 



+ q{u{x, t), v{y, s)) ■ (Va; + ^y)ip{x, y, t, s) 

+ v{A{u{x, t)), A{v{y, s))) gr[^p{-, y, t, s)]{x) 

+ viAHx, t)),A{viy, s))) gri^Pix, •, t, s)](y) 

+ rj{A{uix, t)), A{viy, s))) 5'^[i^(-, ■,t, s)]{x, y) dw > 0. 



Thanks to the regularity of the test function ^, we can now take the limit as r — ^ 
in p.6p . and end up with 



(3.7) 



where 



Qt jjq 



ri{u{x, i), w(y, s)) {dt + ds)ip{x, y, t, s) 

+ q{u{x, t),v{y, s)) ■ (Vj; + Wy)tp{x, y, t, s) 

+ il{A{u{x,t)),A{v{y,s))) g[ij{-,-,t,s)]{x,y) dw>0, 



g[ip{;-)]{x,y)^P.V. 



ip{x + z,y + z) -ip{x,y) 



z >o 



\d+X 



dz. 



Inequality p.7p concludes the first part of the proof. 

2) We now specify the test function ip in order to derive the L^-contraction 
from inequality (13.71) : 



X — V 

1p{x, t, y, S) ^ibpl z I LOn 



i- s 



X + y t + s 



for p > and some (p e C^{Qt) to be chosen later. Here LOp{x) = ujp{xi) ■ ■ ■ ujp{xd) 
and i^p(s) = -^{-) for a nonnegative lo S C^(R) satisfying 



a;(— s) — w(s), uj{s) — for all \s\ > 1, and / uj{s) ds = 1 
The reader can easily check that 

't- s 



{dt+ds)ip{x,y,t,s) ^ujp(^—^ji^p(^—^j{dt + ds)(l>\^- 
(V, + VyMx, y, t, s) = upi^^^^p (^) (V, + Vj,)0( 

9M-,-,t,s)]{x,y) ='^p[-^—)^p[-^-)9[(l}[-,^^) 



2 ' 2 

X + y t + s 



2 ' 2 

x + y 



Note that with this choice of test function V', expressions involving g naturally 
transform into expressions involving g. 
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We now show that, in the hmit p — >■ 0, inequahty l\'3.7\i reduces to 



(3.8) 



ri{u{x, t), v{x, t))dt(t){x, t) 

+ q{u{x, t), v{x, t)) ■ V0(x, t) 

+ r]{A{u{x,t)),A{v{x,t)) g[(l){-,t)]{x) dxdi > 0. 



Loosely speaking the reason for this is that the function ujs converges to the S- 
measure. A proof concerning the local terms can be found in e.g. iSlj. It remains 
to prove that 



M:= 



\A{u{x,t))-A(v{y,s))\ 



'T •'•' WT 



X ~ y \ I t ^ s 

^n I — - — I w„ I ^r— I g 



t + S 



x + y 



dw 



\A{u{x,t)) - A{v{x,t))\ g[(b{;t)]{x) dxdt 
To see this, we add and subtract 



p^O 



\A{u{x,t)) - A{v{x,t))\ 



x — y \ I t — s 
Wo I ^ I ^o I — ^ I g 



}{■■-¥) 



x + y 



dw. 



use the fact that JJ^j ujp {^^) i^p (^) dyds — 1 for any fixed t e (0,r) for p 
small enough, and that </> has compact support in (0,r) to find that 



M < 



Qt jj Q 



QtJJQ 



\A{u{x, t)) - A{v{y, s))\ - \A{u{x, t)) - A{v{x, t))\ 



-p(V)"''(V)-9K'^)l( 

^{■,'-±^)]{^^)-gm;m^) 



dw 



uJo I ^^—T^]uJo ( —r^ 1 \A{u{x,t)) - A{v{x,t))\ dw. 



Let Ml and M2 denote the two integrals on the right hand side of the expression 
above. By the inequality ||a — c| — |6 — c|| < |a — fil we see that 



Mi<K. 



\A{v{x,t))-A{v(y,s))\ Cjp 



x-y 



i- s 



dw. 



since, for aU (a;,i), (y, s) e Qt x Qt, 



(3.9) 



>(.^)] (^) 



<K^:^'±\\D^cj,\\l. 



I<i 



\d+\ 



dz + 2cA||(/'||Lo 



dz 



\z\>l 



\d+\- 
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Note that both mtegrals in l\'3.9\i are finite (use polar coordinates to see this). Using 
the change of variables x — y — h and t — s — t^ we obtain 



M 



\A{v{x, t)) - A{v{x + h,t + T))\ Qjp Q j LOp (0 



<K^ // \Mv(x,t))-A(v(x + h,t + T))\(:b„{-]u}„{-]dxdtdhdT 

<K^II <^p[^)^p{^) [ff \Aiv{x,t))-A{v{x + h,t + T))\dxdt\dhdT 
<K^ sup { ff \A{v{x,t))~A{v{x + h,t + T))\dxdt\^^Q 

\h\,\T\<p\JjQT J 

by continuity of translations in L^ . We refer to Lemma 2.7.2 in [3S] for a similar 
proof. A similar argument using the fact that g[(j)\ G C([0,T]; L^(]R'^)) (cf. Remark 
12. 2p shows that M2 — >■ as p — >■ 0, and we can therefore conclude that M < 
Ml + M2 — > as p — > 0. The proof of p.Sp is now complete. 

3) We now show that inequality p.Sp can be reduced to 
(3.10) // \u{x,t)-v{x,t)\x'{t) dxdt>Q, 



if we take </> = {pr{x)x{t) and send r — > cx) for r > 1, x G ^^(0, T) (with derivative 
x') to be specified later, and 



iprix) ^ u}{x-y)l\y\^r dy. 



All derivatives of ipr are bounded uniformly in r and vanish for all ||a:| — r\ > 1. 
Concerning the fiux-term in (13.81) . we find that 



sgn(u(a;, t) - v{x, t)){f{u{x, t)) - ,f{v{x, t))) ■ \I(j){x, t) dxdt 

Qt 

<Lf\\x\\L^ jj (|u(a;,i)| + |i;(x,t)|)l||,|_,|<idxdr-^0 

by the dominated convergence theorem since u and v belong to L^ and l||^|_,.|<i — > 
as r — >■ c» for all a; G M''. The term in p.Sp containing the non-local operator 
also tends to zero as r — > 00. To see this note that |(7[(^r](a;)| is uniformly bounded 
in r, cf. p.9p . so by integrability of u and v and Holder's inequality, 

\A{u{x,t)) - A(v{x,t))\ \g[vr]{x)\ dxdt 

Qt 

< La(\\u\\li^q^) + \\v\\l1{Qj.)) sup \\g['Pr]\\L--{QT) < '^■ 
\ J r>l 

Hence we find that the integrand is bounded by an L^-function uniformly for r > 1: 

\A(u(x,t)) - A{v(x,t))\\g{tfr\(:x)\ < LA\{u{x,t) ~ v{x,t)\snp\\g[ipr]\\L-°{QT)- 

r>l 

Then for any x, z E M'' fixed and r > |a;| -I- 1, 1^9^(2;) — 1 and 

\ipr{x + z) - (Pr{x)\ < |l|a;+^|<r-l - 1| < l|2|>r-l-|a;| • 

With this in mind we find that 

-*- l^l >r— 1— |a;| j r— foo , 
\z\>0 1^1 

and hence we can conclude by the dominated convergence theorem that 



/[9:'.](x)i< / 'T;;;" dz — 0, 



r-^00 



lim // \A{u{x,t)) - A{v{x,t))\ \g[ipr]ix)\ dxdt = 0. 



12 S. CIFANI AND E.R. JAKOBSEN 

4) To conclude the proof, we now take X — Xm ^^^ 

Xp(i) = / (a;p(T - h) - uj^{t - is)) dr, 

^ — oo 

where r > 1 and < ti < t2 < T. Loosely speakmg, the function X/i is a smooth 
approximation of the indicator function l^f^ j,) which is zero near t = and t = T 
when /z > is small enough. Since x^(^) ~ ^fi{t — ti) — ^t^it ~ ^2), inequality p.lOp 
reduces to 

|M(a;,i) — v(a;, t)| a;p(i — ^2) dxdi < // \u{x,t) — v{x,t)\ujfj,{t — ti) dxdt. 

jJqt 

By taking /i small enough and using Fubini's theorem, we can rewrite this inequality 
as 

(3.11) <i>*uj^,{t2)<^*ujf,{ti) for ^{t}^ f \u{x,t)-v{x,t)\dx, 

where (/>i * 02(t) = J^(f>i{s) (f>2{t - s)ds. Since u,v e C([0,T]; i^(M'')), we see that 
$ e C([0,r]), and hence by standard properties of convolutions, 

$ * Wp(t) -> $(t) as ^ -> 0. 

for all t G (0, T). Hence we can send /i — > in p. lip to obtain 

||(u-i;)(-,i2)||Li(R<i) < ||(w- w)(-,ti)||i,i(Rd). 

Finally, the theorem follows from renaming i2 and sending ii — >■ using Hi) and 
C{[0,T];L^{R'^)) regularity of u and v. D 

4. A CONVERGENT NUMERICAL METHOD 

In this section we introduce a numerical method for the initial value problem 
(II. ip which is monotone and conservative. Then we prove that the limit of any 
convergent sequence of solutions of the method (as Ax -^ 0) is an entropy solution 
of (jl.ll) . Finally we prove that any sequence of solutions of the method is relatively 
compact whenever the initial datum is a bounded integrable function of bounded 
variation, and hence we establish the existence of an entropy solution of (II. ip in 
this case. Some numerical simulations based on this method are presented in the 
last section. 

4.1. Definition and properties of the numerical method. For simplicity we 
only consider uniform space/time grids and we start by the one dimensional case. 
The spatial grid then consists of the points Xi = iAx ioi i ^ 1 and the temporal 
grid of tn = nAt for n = 0, . . . ,N and NAt — T. The explicit numerical method 
we consider then takes the form 

jjn+i ^jjn_ ^tD-F{U^, f/f+i) + AtY,G,{AiU-+,) - A{Un), 
U-^^ I Mx) dx, 

^X Jxi+Ax[0,l) 

where D~Ui = -^{Ui — Ui-i), F : M^ ^- M is a numerical flux satisfying 

a) F is Lipschitz continuous with Lipschitz constant Lp, 

b) F is consistent, F{u,u) — f{u) for all w S M, 

c) F{ui,U2) is non-decreasing w.r.t. ui and non-increasing w.r.t. U2, 
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and Gi is defined by 

Gi^ cx I I .... for i ^ 0. 

In tlie mufti dimensionaf case tfie spatiaf grid is Aa; J/^ (Ax > 0) witir points 

Xa = Ax a where a = (ai, . . . , ad) G Z . 

Let ei be the d-vector with /-component 1 and the other components and define 
the two box domains 

R^l^x[Q,lY and Ro = ^[~\,lf, 

noting that Ua{xa + R) — ^a{xa + Rq) — Hi''- The exphcit numericai method we 
consider now tafces the form 



(4.1) <^ 



d 

[/^i = U-- AtY^D^Fiiu:, U:^J + AtY,Gp{A{u:^^) - A{u:)), 
1=1 p^o 

^°=A^ "o(a;) dx, 



where D^ Ua = -^{Ua — Ua-ei), -F) : M^ — > M is a numericai flux satisfying a) - c) 
above with /; replacing /, and Ga is defined by 

f dz 

Ga=cx g^j^ for a^O. 

Jx^ + Ro 1^1 

Note that Ga is positive and finite since ^ Xa + Ra unless a = 0. 
Remark 4.1. An admissible numerical fiux Fi is e.g. the Lax-Friedrichs flux, 

We refer the reader to f23] or |26[ Chapter 3] for a detailed presentation of more 
numerical fluxes which fulfill assumptions a) - c). 

Let us introduce the piecewise constant space/time interpolation 

u{x,t) = Ua for all {x,t) e {xa + i?) x [t„,i„+i). 

In the following we often need the relation 

/.^N \^r (A(TT^ \ A(TT^\\ n f A{u{ya+Z,tn))~A{u{ya,tn)) ^_ 

(4.2 2^(^l3(A[Ua+fs)-A{Ua)) = cx r-^-j^ dz, 

where ya = Xa -\ — ^(1, ■ • • , !)• Note that this is an approximation of the principal 
value of the integral since Rq ^ as Ax — ^ in a symmetric way. 

We now check that the numerical method (|4.1I) is conservative and monotone. 

Lemma 4.2. The numerical method (14.11) is conservative, i.e. 
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Proof. First wc show that Y.aez" |C/q | < oo for all n = 0, . . . , TV. By dUT]), 
(4.3) 

E \u:+'\< E {\u:\ + Atj2\DrFi{u:,u:^j\ 

^ E { i^"i + I^E (^^ |t^a - c/^-j + ^i=^ K, 

aeZ'' ^ 1 = 1 

+ Ai^G^(|^(C/„\^)| + |A(C/:)f 



- U"\ 



'■ dz 



^#0 
where, using that {z : \z\ < ^} C Rq, 

Since Ax"* X^aez^' l^al — li'"o||Li(R<i) < oo, we can iterate estimate (|4.3p to find that 
X^aGzrf |C^a I < oo ^nd hcncc lini|Q|^oo \Ua\ = for all n = 0, . . . , iV. 
Now we sum (14. ip over a to find that 

d 

E ^^'= E c/:-At^ ^i?r^^(c^a",c^:+e,) 

aez<* QeZ"* aei'' i=i 

aeZ<* ^7^0 

The proof is now complete if we can show that the F and G sums are equal to zero. 
The F-sum is telescoping and since lim|(j|_j.oc |C/^| = 0, 

E D7F,iu:,u-^j = E ^^'""'''"^"^^/^''"'^"'""^ = 0. 

To treat the G-sum, note that we have found above that 

^ ^ G^ \A{u:^f,) - A{U:)\ < 2LAAtY^Gp ^ \U2\ < oo, 

a6Z<i P^O P=jt0 QSZ'i 

and we also have that X^a I^(^q)I — FA^aWa\ < oo. In view of this we can now 
change the order of summation, and split the sums to find that 

Y, ^ G^{A{U-^^) - A{U-)) = E G^ E {^U'^+p) - MU^) 

= E G/3 ( E ^(f^^) - E ^(^"")) = 0- 

The proof is now complete. D 

Next, we check monotonicity by showing that the right-hand side of the numerical 
method (I4.ip is a non-decreasing function of all its variables Ua. This is clear for 
all Ua such that (3 ^ a since the numerical flux Fi is increasing w.r.t. its first 
variable, non-increasing w.r.t. its second one, the function A is non-decreasing, and 



DEGENERATE FRACTIONAL ORDER CONVECTION-DIFFUSION EQUATIONS 15 

the weights G^ are all positive. Then we differentiate the right hand side of (I4.ip 
w.r.t. U^ and find that it is non-negative provided the following the CFL condition 
holds, 

Ai / , f dz \ At 

(4.4) 2dLF— + [c,2'LaJ ^^)^^ < 1. 



I>i i.K+^/Aa;^ 



We have thus proved the following result: 



Lemma 4.3. The numerical method (j4.ip is monotone provided that the CFL 
condition (14.41) is assumed to hold. 



In what follows, the CFL condition (|4.4p is always assumed to hold, and mono- 
tonicity is thus always ensured. 

4.2. Convergence toward the entropy solution. We prove that any limit of 
a uniformly bounded sequence of solutions of the numerical method (|4.ip is an 
entropy solution of (|1.1[) . 

Theorem 4.4. // {u} is a sequence of solutions of (j4.ip . uniformly bounded in 
L^{Qt), and there exists u G L°°{Qt) n C([0,T]; ^^(M'*)) such that u ->■ u in 
C{[0,T]; L^{W^)) as Ax — >■ 0, then u is an entropy solution of (|l.ip . 

Proof. Note that part i) in the definition of entropy solution (Definition 12. ip is 
already satisfied. Part Hi) follows since ||u(-, 0) — uo||ii(Rd) -> as Ax — > by the 
definition of u. What remains to prove is part ii). 

First we prove that the numerical method (|4.ip satisfies a discrete entropy in- 
equality which resembles the one in ii), Definition [5TTJ To this end, let us introduce 
the notation a Ab = min{a, b} and a\/ b = max{a, 6}, choose an r > 0, and exploit 
monotonicity to obtain the inequalities 

d 

C/^+i Vk<U2Vk-AtJ2 DiFi{U2 V fc, U2+e, V k) 

1=1 

+ At J2 ^P {MK+p V fc) - AiU^ V k)] 

0<Ax\l3\<r 



At 1(,,+^) ([/„"+!) Yl Gp(^A{U:^p)-A{u:) 



Ax\P\>r 

and 

d 
[/:+! Ak<U^Ak-AtY, DiFi{U2 A k, U2+,, A k) 

1=1 

+ At Y. ^P {a{U^+P a A) - A{Ua A k)] 

0<Ax|/3|<r 

+ Atl(_^,,)(L/ri) Y. G0(a{U-^^)~A{U^^ 

Ax\P\>r 

Note that the polygonal set 

Pr := U {xp + Ro) 

0<Ax\p\<r 

{xj3 = Ax 13) does not include points from the box i?o, and converges to the punc- 
tured ball {2: : < |z| < r} as Ax — ?► in the sense that lp,X-z) — > '^Q<\z\<r{z) a.e. 
as Ax — )■ 0. 

Remember that rjkiUJ^) = \U2 — k\, and let 

QhAK) = Fi{u2 V fc, U2+,, V fc) - Fi{U2 A fc, [/:+,, A k). 
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Thanks to the relations 



|u — fc| — u\/ k — u A k, 
\A{u) - A{k)\ = A{u V fc) - A{u A k), 



we can subtract the above two inequahties to obtain that 



1=1 



0<Ax\l3\<r 

Atry;(t/:+i) Yl G0[a{U:+p)-A{U:))<O 

Ax|/3|>r 



Let us take a nonnegative function ip € C^{Qt), and define ip"^ ~ (p(xa,i„). If 
we multiply both sides of the above inequality by ip2 ) sum over all a G Z'* and all 
n G {0, . . . , N}, and use summation by parts for the local terms, we end up with 
the cell entropy inequality 



N d n n 

+ Ax'AtY: E Eom(^:) ^"^" ' ^° 



(4.5) 



Ax 

n=OaeZ'' 0<ArE|^|<r 

A:r^AiEE^^(f^^')^S E G,[a{u:^^)-A{u:))>0. 

n=0 ai^ld A2;|;3|>r 



To derive this inequality we have used the change of indices (/3, a) ~> (— /?, a + (3) 
to see that 



AT 

Ax'^AtE E ^S E G^(^A(fc)(A(C/„\^)) -ry^(fc)(A(C/:))) 

n=OQeZ'' 0<A£i;|/3|<r 

AT 

= Ax"^ Ai E E ^^Mk) (MUn) E G^ (^S+^ - ^S' 

«=OQgZd 0<A2;|;3|<r 
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Let Ra — Xa + R- We now claim that for each fixed Aa; > 0, inequality (I4.5P 
implies 



(4.6) 

n=0 aGZd 1=1 ''*" "^ ^^ ^ 

+ EE'?fc(^^') E G;5(A(f/:+,)-A([/:)) / / ^{x,t)Axdt>Q. 



To see this we proceed by contradiction, and assume that (|4.6p is strictly negative. 
We then sum together several inequalities of the form (|4.5p where, instead of (p^ = 
f{xa,tn) which are computed on the original space/time grid (xQ,,i„), we use the 
values ip^ = ip{xa,t„) computed on the finer grid {xa,tn) where Xa = {Ax/M)a 
while tn = n{At/M) for some M > 0. Note that, since all these inequalities of 
the form (j4.5p share the same underlying numerical solution (C/"), they can be 
rearranged as one inequality, i.e. 




Vj - f^ 



in — l 

M ) M ^ ^ 



EEEQ^ra (t)'^ Y. Y. 

n=OaeZ'^l=l \^ ' m:t„e[t„,t„ + i)7:«",6-R. 




E E '?^w(^(f^")) 

Aa;\'' Ai 



f ) f E E E ^.(^."V.-v'T 

a 



"=OqgZ'' Ax|/3|>r 

m: tme[t„,t„+i) 7:a;Teficv / 

(loosely speaking, by summing all these inequalities of the form (j4.5p together we are 
filling the mesh-sets Ra x [i„, t„-|-i) with several samples of the test function (p\ this 
has been done in order to recreate in each mesh-set a Riemann sum approximation 
which gets closer and closer to its respective integral as the value of the control 
parameter M increases). The Riemann sum approximations in the first, second, 
and fourth term of (|4.7p are arbitrarily close to their respective terms in (|4.6p as 
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M increases. For the third term in (|4.7p note that, cf. 

^ m:t„e[t„,t„+i)7:S7e-Ra 0<Aa;|^|<r 

(4-8) , d ' . . 

' Ax\ ^ Sr^ \^ f 'ig(;/7 + Z, tm) - IfiiVf, tm) , 

Ml M ^ ^ J,^p |zK+^ 

(the definitions of (^,2/7 are analogous to those of u,ya) and so the Riemann sum 
approximation on the right-hand side of (J4.8I) is, as M increases, arbitrarily close 
to 



/ rn^x dzda;dt. 



|^|d+A 



This is due to the fact that, since we are integrating away from the singularity, the 
right-hand side of (|4.8p is well defined, and the sum over all {x-y, tm) can be moved 
inside the integral z ^ Pr. Therefore, since (I4.7P is arbitrarily close to the left-hand 
side of (|4.6p . the left-hand side of (|4.6I) cannot be negative, and we have produced 
a contradiction. 

Using the piecewise constant space/time interpolation u, we can now rewrite 
inequality (|4.6p as 

< rik{u{x,t))dtip{x,t) +y2Qh.i{u{^,t)) da;i(p{x,t) 

fAf-f ^\\\ f V{x + Z,t)-ip{x,t) 

(4.9) + ^^(^-) (^("(^' ^))) ^ ^pTA dz 

,, , A NN / X /" A(u(a; + z, t)) — A(u(x, t)) , , , 

+ f]i{u{x,t + At))^(x,t) -^-^ '/ "- ^ ' '' dz\ dxdt 

Jr''\p^ \z\ ^ J 

>0{Ax)+0{At). 

Convergence up to a subsequence for the first three terms in (j4.9p is immediate 
thanks to the a.e. convergence of u toward u. For the local terms this is already 
well known, cf. |26l Theorem 3.9]. For the term containing the inner integral Pr, 
convergence follows thanks to the convergence of Ip^ — > lo<|z|<r a-e., the proper- 
ties of ip {Jp '^ u\d+\ dz is uniformly bounded and compactly supported), 
uniform boundedness of u, and the fact that the function r]k{-) is continuous. 

To conclude, we need to establish convergence for the term containing the dis- 
continuous sign function ri'f.{-), and we argue as in [19] (p. 109). First note that 
since u ^ u in C{[0,T];L\m'^)), also u{-,- + At) ->■ u in C{[0,T];L\R'^)) and 
a.e. for a subsequence. Then note that 77^ (s) is continuous for s ^ k, and that the 
measure of the set 

Uk^{{x,t) gQt: u{x,t)^k} 

is for a.e. fc G R. For such k, ?7^(u(- 4- At, •)) -^ v'ki'^) ^-e., and we can go to the 
limit in the term involving 77^ in (|4.9p using the dominated convergence theorem, 
177^1 ^ Ij and uniform boundedness of w and A(u). 

For the remaining k, we use an approximating sequence made of those k for 
which convergence holds true. To be more precise, let am, hm be sequence of values 
such that meas(Z//a„i) = nreas(Z-/b^) = 0, where «„ /* k and bm\i k. Note that the 
mean value 

2('?L('") + '7L(")) ^V'kiu) asam,bm ^ k. 
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Thus we can use the entropy inequahty for the sequence Um and the entropy in- 
equahty for the sequence bm, take the average, and go to the hmit to prove the 
entropy inequahty for every critical value k. Convergence for the whole sequence u 
is a consequence of uniqueness for entropy solutions of (jl.ip . D 

4.3. BV initial data: Compactness and existence. We now show that the 
sequence of solutions of the method, {u : Ax > 0}, is relatively compact whenever 



Using this result and Theorem 14.41 we then obtain existence of an entropy solution 
of the initial value problem (|l.ll) . We start by the following a priori estimates. 

Lemma 4.5. If uq G L°°{R'^) D L^{R'^) D BV{R'^), then, for all t,s>0, 

i) II"(-iOIIl~(R'') < ll"o||L°=(R<i)> 
i^) II"(-)0IIl1(R'') < ll"o||Li(R'i), 

in) |w(-,i)|BV(R'i) < |wo|Bv(Rd), 

iv) \\u(-, s) — ■u(-, t)llii(Rd) < a(\s — ^1 + Ai) where, for some c> 0, 

{c\s\ */Ag(0,1), 

cr(s) = < c|slns| if X — 1, 

[c\s\i */Ag(1,2). 

Lemma 14.51 along with a Kolmogorov type of compactness theorem, cf. The- 
orem A. 8 in [55], yields the existence of a subsequence {u} which converges in 
C([0, T]; Ll^^{R'^)) (and hence a.e. up to a further subsequence) toward a limit u as 
Ax — > 0. Moreover, the limit u inherits all the a priori estimates i)-iv) in Lemma 
14.51 (with Ai = 0). Moreover, by ii) and the dominated convergence theorem, we 
see that u ^ u also in C([0, T]; L^(M'')). In short, we have the following result: 

Lemma 4.6. The numerical solutions {u : Ax > 0} converge, up to a subsequence, 
toward a limit u in C{[0,T]; L^{R'^)) as Ax — >■ 0. Moreover, 

u e L°°{Qt) n C{[0,T];L\m.'^)) n L°°(Q,T; BV{m.'^)). 

Lemma 14.61 and Theorem 14.41 implv the following existence result: 

Theorem 4.7. (Existence) If uq G L°° {W^) D L^ {R'^) , then there exists an entropy 
solution of the initial value problem (jl.ll) . 

Proof For initial data uq G L°°{W') n L^{R'') n BV{R'^) existence is granted by 
the numerical method (|4.ip (Lemma 14. 6p . For more general initial data uq G 

L°°(R'*) n Li(M''), we consider approximations uo,„ G L°°(M'^) n ^^(M'') n BV{R'^) 
such that 

11^0 - uo,n||Li(Rd) ^- as n -^ oo. 

Let Um,Un denote the entropy solutions corresponding to uo,n)Wo,m respectively, 
and use the L^-contraction (Theorem 13. ip to see that 

\\Un - '"m||c([0,T];Li(R'')) < li"0,n " ""Cm |i Li(R<i) -> aS 71, 771 -> OO. 

Therefore, the sequence of entropy solutions {w„} is Cauchy in C([0,T]; L^(R'^)) 
and admits a limit u. To prove that u is also an entropy solution of (|l.ip . one can 
pass to the limit 77 — ^ oo in the entropy inequality for Un- D 

Proof of Lemma \4.5\ The maximum principle i) is a direct consequence of mono- 
tonicity. To see this let s — sup^g^ti Wa\^ and choose C/" = s to obtain that 

d 

[/„"+! < s - AtJ^D^Fiis, s) + AtJ2 G0{A{s) - A{s)) = s. 

1=1 13^0 
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Similarly, choosing C/" = — s, one obtains U^~^^ > — s. Furthermore, since the nu- 
merical method (j4.ip is conservative, monotone, and translation invariant (trans- 
lation invariance is a consequence of the fact that the numerical method does not 
explicitly depend on the variables Xa,tn), inequalities ii)-iii) are consequences of 
the results due to Crandall- Tartar [ TS] . 

We now prove iv). By (14.11) and Lipschitz continuity of Fi, 

1=1 ^ ^ p=^o 

Let us multiply by Ax'' in the above inequality, and sum over a € Z'^ to see that 

Ax'' Y, |c/:+'-c/:| 

d 

< 2LpAx''~^AtJ2 E l^"+e. - U2\ + Ax" At ^ ^ Gp\A{U-+^) A{U2)\. 
Let «"(•) = u{-,tn) and note that the first term then is equal to 

d „ 

2LFAty] / |u"(-,x')|sy,,(M) dx' < 2dLpAt\u"\BviR-i) - 0{At), 

while the second term can be estimated by (cf. (|4.2p ) 

Ax''Y.llGAAiU2+p)~A{U2)\ 

\J^<\z\<l J\z\>\J ^^^^ Fl 

Easy computations in polar coordinates show that the second integral is 0(1) while 

. I I fO(l) if AG (0,1), 

/a.= / J|i3^dz=<^0(|lnAx|) ifA-1, 

•^^<l^l<^'' \0{Ax^-^) if A e (1,2). 

Adding all the above estimates yields 

Ax'' Y. |[/:+i - Ul\ = 0{At) + 0{AtlA.) + 0{At). 

By the CFL condition (|i^ . A^/aj; = cr(At), and the resuh follows. D 

5. Extension to general Levy operators 

The ideas developed in this paper can also be used to establish well-posedness 
for entropy solutions of a more general class of fractional equations of the form 

{dtu + V-f{u)=g^[A(u)] inQT-K''x(0,T), 
1u(x,0) =U(i{x) \tlW\ 
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where the fractional Laplacian g has been replaced with a more general Levy op- 
erator g^: 



9^^[4>]{x) = / 0(x + z) - <l){x) - z ■ V0(x)l|;,|<i d^(z), 

^|z|>0 

where the Levy measure /i is a positive Radon measure satisfying 



(5.2) / |z|^Al ^(dz) < oo. 

-'|z|>0 

Note that g^ is self-adjoint if and only if /i is symmetric: n{—B) = ijl{B) for all 
open sets B. The adjoint 5* (defined through ^ ug^[v\ = ^ g*p[^]v) equals 

9lW\{x) = / (j){x- z)- (j){x) + z ■ V<?!)(a;)l|2|<i d^(z). 

A Taylor expansion shows that both 5^ [0] and g* [(j)] are well defined whenever (j) is 
C^ and bounded. The gradient term is needed when /i is not radially symmetric, in 
the radially symmetric case g^, (= g*) can be defined as before as a principal value 
and no gradient term. The operator g^ is the generator of a pure jump Levy process 
and these processes have many applications in Physics and Finance, cf. e.g. ^. 

We need a modified definition of Entropy solutions. Remember the notation rjk 
and Qk introduced in Section [2j and define for r > 0, 

where g^^rM{x) = g^[Lp{-)l\^\<^]{x), 



9lM{x)=l if{x + z)-Lp{x) ^.{Az), and l1=l z ^i{dz). 

J \z\ >r J r< l^l <1 

We also use the notation g* ^ and ^!^'* for the adjoint operators, and note that 

Let us point out that the adjoint operator g* could have also been defined as g^ 
with i'{B) = fi{—B). From this equivalent definition it is clear that the adjoint 
operator g* is still a Levy operator. 

Definition 5.1. A function u is an entropy solution of the initial value problem 
(|5.ip provided that 

i) ueL°"(QT)nC([o,T];Li(R'')); 

ii) for all fc € R, all r > 0, and all nonnegative test functions Lp € C^{Qt), 

f]k{u)dt(p + qk{u) ■ Vip + ■qA(k){A{u)) 9l rM 
Qt 

+ V'kiu) g''A[A{u)] ^ + TiA(k){A{u)) 7; • V(p dxdi > 0; 

Hi) u{x,(S) — uq{x) a.e. 

Remark 5.1. All terms in ii) are well-defined in view oii). Except for the gj^-term, 
this follows from the discussion proceeding Definition 12.11 - see Remark 12.21 Note 
that the integrand oi g^^[A{u)] is measurable w.r.t. the product measure d/i(z)dxdt 
since since it is the d/x(z)da;di-a.e. limit of continuous functions. This follows readily 
from the fact that u is the dccdi-a.e. limit of smooth functions. Integrability then 
follows by Fubini's theorem, integrate first w.r.t. to dxdi and then w.r.t. diJ,{z) 
using dOl)- By Fubini we also see that glXA{u)] G C([0, T]; ^^(R'^)) and it easily 
follows that the gj^-tcrm is well-defined. 
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Again classical solutions are entropy solutions and entropy solutions are weak 
solutions. The proof is essentially the same as the one given in Section [2] with the 
additional information that whenever A{u) is smooth 

Vk{u{x))S/[A{u{x))] ^ r,'^^,,){A{u{x))V[A{u{x))] ^S/[rjAikM{u{xm a.e. 

We also have a i-'^-contraction and hence uniqueness result: 

Theorem 5.2. Let u and v be two entropy solutions of the initial value problem 
(IS.ip with initial data uq and vq. Then, for all t G (0,r), 

ll"(-:*) -^v{-,t)\\m^d) < ||uo-uo||li(R'')- 

Proof. We proceed as in the proof of Theorem l3.1l let us take the entropy inequality 
for u = u{x,t) and the one for v = v{y,s), integrate both in space/time, and sum 
the resulting inequalities together to obtain an expression equivalent to p.3p . At 
this point we use the change of variables {x,y) — ?> {x — z,y ~ z) to obtain the 
inequality 

r]{u{x, t),v{y, s)) (dt + ds)'ip{x, y, t, s) 

+ q{u{x, t),v{y, s)) ■ (Vj, + \^y)'ip{x, y, t, s) 

+ 7]{A{u{x, t)),A{v{y, s))) g;^M-,y, t, s)]{x) 

+ T]iA{u{x, t)),A{v{y, s))) g;j^{x, •, t, s)]{y) 

+ r]{A{uix, t)),Aiviy, s))) 7^ • (V, + Vy)ipix, y, t, s) 

+ r]{Aiuix,t)),A{viy,s))) gl;*[i>i-,-,t,s)]{x,y) dw>0, 

where 

9l^*M-, ■)]{x, y) = / (fix- z,y- z)- Lp(x, y) ^i{dz). 

J\z\>r 

We can now send r — >■ and recover the equivalent of expression p.7p in the present 
setting, 

■q{u{x, t), v{y, s)) {dt + ds)ip{x, y, t, s) 
IQt JJQt 

+ q{u{x, t),v{y, s)) ■ (Vx + Vj^)V'(a;, y, t, s) 

+ r,iA{u{x, t)),A{v{y, s))) g;[i^{; •, t, s)]{x, y) Aw > 0, 

where 

9tl'Pi-r)]{x,y) 

= / 'fiix - z,y~ z)- tp{x, y) + z- (V^; + Vy)i^(a;, y)l|^|<i ^^(dz). 

"'|2|>0 

From now on, the proof follows the one of Theorem 13.11 (just replace the operator 
g therein with the operator g*). D 

Existence of solutions can be obtained e.g. by the vanishing viscosity method 
and a compensated compactness argument, but we do not give the details here. We 
just remark that the vanishing viscosity equations have smooth solutions since the 
principle term is the (linear 2nd order) Laplace term. 

Theorem 5.3. There exists a unique entropy solution of the initial value problem 

(EH). 
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6. Connections to HJB equations 

In one space dimension it is well known that the gradient of the (viscosity) 
solution of a HJB equation is the (entropy) solution of a conservation law, see 
e.g. [3S]. Variants of this result are still true in the current fractional setting as we 
will explain now. First we consider the following two initial value problems in one 
space dimension: 

j ut + f (dxu) + g[u] ^ edlu in Qt, 
1 u{x,0)—uo{x) in R, 



(HJB) 

and 

(FCL) 



vt + dxf{v) + g[v\ = edlv in Qt, 

v{x, 0) = dxUo{x) in M, 

for any e > 0. The first equation is a HJB equation and the second one a frac- 
tional conservation law. To simplify, let us consider the following strong but rather 
standard regularity assumptions: 

(al) / G C2(R), 

(a2) uq is bounded and Lipschitz continuous, and 
(a3) d^uo is bounded and belongs to ^^(R) n BV{R). 
Standard results then show that: 

(i) there is a unique bounded Holder continuous (viscosity) solution u^ of 

(HJB) foranye>0 [H], 
(ii) there is a unique bounded (entropy) solution v'^ G i^(0, T; BV n L^) of the 

(fractional) conservation law for any e > [U [H] , 
(iii) when e > both u^ and v^ are C^, 
(iv) u^ —> u° uniformly [35] and w^ —> v'^ in L^ [T] as £ -^ 0. 

By differentiating (HJB) and using uniqueness for (FCL), we find that 

for any e > 0, and hence 

v^(j) ^ - ff u'd^ct) for any G C^iQr)- 



Sending e — ^ in the above inequality using dominated convergence theorem and 
(iv) then leads to 

^^V = -^^i"a,'^ for any cj^eC^iQr), 

and we have the following result: 

Theorem 6.1. The distributional x-derivative of the viscosity solution of (HJB) 
is equal to the unique entropy solution of (FCL). 

The only part missing in the proof of this theorem, is the proof of (iii). This 
result follows e.g. from energy estimates and standard parabolic compactness re- 
sults (yields L^{0,T;H^) solutions) combined with regularity theory for the Heat 
equation, interpolation, and bootstrapping arguments (yields smooth solutions). 
We skip the long and fairly standard details. 

If we drop the convection term, we get a similar correspondence in any space 
dimension. Consider the following two initial value problems: 

( ut- A{gf_,[u]) ^ eAu in Qt, 
(HJB2) <^ 

I u{x,0)=uoix) inR^, 
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and 

f vt - 5p[A(w)] = eAw in Qt, 

(FDE) <^ , 

1 w(x,0)=.g,jMo](a;) inM^ 



for any e > 0. The first equation is still a HJB equation while the second one is a 
degenerate fractional diffusion equation. To simplify, let us consider the following 
rather strong regularity assumptions: 

(bl) A S C^(R) is non-decreasing and Lipschitz continuous, 
(b2) lio is bounded and Lipschitz continuous, and 
(b3) g^[uo] is bounded, BV, and belongs to L^ . 
Again we have the following of properties: 

(i) there is a unique bounded Holder continuous (viscosity) solution u^ of 

(HJB2) for any e > 0, 
(ii) there is a unique bounded (entropy) solution u^ e L^{Q,T\BV C\ L^) of 

(FDE) for any e > 0, 
(iii) when e > both m^ and f^ are C^, 
(iv) u^ -^ u'^ uniformly and u^ ^- u° in L^ as e — >■ 0. 

By applying 5^ to (HJB2) and using uniqueness for (FDE), we find that 

for any e > 0, and hence 

jjv'cj,^jlu'g;[^] for any </. G ^^(Qt). 

Sending e — > in the above inequality using the dominated convergence theorem 
and (iv) then leads to 

^V = ^A;M for any c^eC^iQT), 
and we have the following result: 

Theorem 6.2. // u is the unique viscosity solution of (HJB2), then v = 5/j[w] 
(where g^ is taken in the sense of distributions) is the unique entropy solution of 
(FDE). 

Proof. For the HJB equation well-posedness of viscosity solutions for e > and the 
uniform convergence u"^ — > u'^ is fairly standard and can be found e.g. as a simple 
special case of results in [5S] . 

Existence and uniqueness in (ii) follow from this paper for e = 0. The arguments 
in this paper can easily be extended to include the eAw-term (this is standard) and 
hence we have (ii) for any e > 0. The limit v"^ — >■ v can be obtained through 
a standard Kuznetzov type argument, cf. [H [10] for the case when A is linear. 
We will give the result for the non-linear case in a future paper. The regularity 
for e > is clear since the eAw-term is the principal term in the equation. It 
follows e.g. from (i) energy estimates and a classical parabolic compactness ar- 
gument (yields L^(0,T; i/^(R'^))-solutions) and (ii) regularity theory for the Heat 
equation combined with bootstrapping (yields smooth solutions). The fractional 
term is always related to integer order derivatives through interpolation estimates. 
The detailed proof is long and rather classical and is best left to the interested 
reader. D 

Remark 6.3. Such correspondences between HJB equations and degenerate convec- 
tion diffusion equation can be useful for at least two reasons. 
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1) They allow for integral representation formulas for the solutions of the 
degenerate convection diffusion equations via representation formulas for 
the solutions of the HJB equations. See e.g. chapter 3.4 in ^24^ for the case 
of one dimensional scalar conservation laws. 

2) They allow for efficient numerical methods for the non-divergence form 
HJB equation, by solving the divergence form degenerate convection dif- 
fusion equation by finite elements or spectral methods and then using the 
correspondence (and the HJB equation) to find the HJB solution. 

The solutions of the above HJB equations are value functions of suitably defined 
stochastic differential games (see e.g. |^), i.e. they have integral representation 
formulas. Since HJB equations are fully non-linear non-divergence form equations, 
it is not natural or easy to solve them directly by well-established, flexible, and effi- 
cient methods like the finite element and spectral methods. Such methods do apply 
to divergence form equations like the degenerate convection diffusion equations (cf. 

e.g. [13 1301 HI])- 

7. Numerical experiments 

We conclude this paper by presenting some experimental results obtained using 
the numerical method (|4.ip with d = 1. We simulate fractional strongly degenerate 
equations and compare them to fractional conservation laws and local convection 
diffusion equations. Our simulations give some insight into how the solutions of 
these new equations behave. Note that this type of fractional equations have never 
been simulated (or analyzed) before. 

In our computations, we restrict ourselves to the bounded region ft — {x : \x\ < 
2} and impose zero Dirichlet boundary conditions on the whole exterior domain 
{x : \x\ > 2}. We consider the degenerate fractional convection-diffusion equations 
with Burgers type convection (/(u) — u^/2), 

(7.1) dtu + udxU ^ g[A{u)], 
and fractional degenerate diffusion equations (/ = 0), 

(7.2) dtU^g[A{u)], 

for two different strongly degenerate diffusions, defined through two different A's: 

Ai{u) — max(u,0) 

and 

r u<0.5, 

A2{u) = < 5(2.5ii-1.25)(u-0.5) 0.5 < m < 0.6, 
[ 1.25 + 2.5(u- 0.6) u> 0.6. 

The numerical experiments below show e.g. how solutions of (jl.ip can develop 
shock discontinuities in finite time for all A e (0, 2). Furthermore, they show that, 
contrary to the linear case, equation ()7.2p does not have smooth solutions for i > 
when the initial data is non-smooth. We also observe that for A ~ 2, solutions are 
very close to solutions of the corresponding local problem with A = 2. 

In figure Figure [T] (a)-(b) we plot the solutions of (|7.ip with linear and non- 
linear fractional diffusion {A{u) = u and A = Ai) to show how the non-linearity 
influences both the shock size and speed. 

Figure [5] (a) shows that a shock discontinuity develops in finite time in the 
region where A2 is zero. This phenomenon is well known for degenerate convection- 
diffusion equations (|1.2p as shown in Figure [5] (b). Here and in what follows, we 
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Figure 1. Numerical solutions of (|7.ip at T ~ 0.5 with Ax = 
1/500 and piecewise constant initial data: (a) A = Ai (solid) and 
A{u) = u (dotted) for A = 0.5; (b) A = Ai with A = 0.001 (dotted), 
A = 0.5 (dashed), and A = 0.999 (solid). 



have used the convergent numerical scheme (cf. [25 ) 



(7.3) U, 



n+l 



C/f - AtD~F{U^, C/f+i) + {2TTfAtD- 



A{u:xi) " ^(un 



'+'" ' ^"■■' ^ Ax 

to compute the solutions of degenerate convection-diffusion equations (jl.21) , 




(a) 



(b) 



Figure 2. Burgers's flux and A = A2 with Ax = 1/500 and 
piecewise linear initial data: (a) solutions of (14.11) with A — 0.3 at 
T = 0.25 (dotted) and T = 0.5 (sohd); (b) solutions of ([731) at 
T = 0.01 (dotted) and T = 0.025 (solid). 

Figure [3] (a) displays the solutions of (|7.2p with A{u) = u and A = A^. Note 
that, when A — A2, the initially discontinuous solution becomes continuous in finite 
time but not differentiable. In the non-degenerate case, dtu = g[u\, the initially 
discontinuous solution becomes smooth immediately for all values of A, cf. Figure 
[3](b). This behavior agrees with results from [5D]. 

In Figure m we compare the solutions of (|7.2p for A « 2, with the solutions of 
a properly scaled equation (|1.2p (A = 2). We use our scheme (|4.ip to compute the 
first set of solutions, while scheme (|7.3p is used to compute the second. Again, we 
have restricted our computational domain to fi. As expected, the solutions of the 
two equations are very close since —(—A) 2 — >• A(j) as A — )■ 2 for regular enough 
(j). The two methods are however fundamentally different: (|7.3I) uses a three-points 
stencil, while (|4.1I) uses a whole-domain stencil. 
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Figure 3. Solutions of ([73) with Ax = 1/500, A = 0.3, and 
piecewise constant initial data: (a) A — Ai with T = 0.1 (dotted) 
and T = 3 (sohd); (b) A{u) = u with T == 0.1 (dotted) and T = 3 
(sohd). 





Figure 4. A{u) = A^ with T = 0.005, Aa; = 1/500, and piecewise 
constant initial data: (a)-(b) solutions of the non-local numerical 
method (|4.ip (solid) with A w 2 compared with solutions of the 
local numerical method (|7.3p (dotted). 



Appendix A. A technical result 
In this section, we prove a technical result used in the proof of Lemma 
Lemma A.l. Let u G BV{R'^), then 

(A.l) / \u{x + z) ~ u{x)\ dx < vd\z\\u\Bv(M'')- 

Note that a more refined argument would give a factor 1 instead of \/d in (lA.ll) . 
This is unimportant in this paper and we skip it. We now give a proof for (lA.ll) in 
the case d — 2, analogous ideas can then be used in higher dimensions. 

Proof. We define the total variation |u|By(R2) as, cf. [2B1 expression A. 19], 

(A.2) |u|By(R2) = / \u{xi,-)\BviR) (1x1+ / I"(->2^2)|bu(r) da;2. 

Jm JR 

Then, since /jj \u{x + z) — u{x)\ dx < \z\\u\gY^g^), we write 



\u(x + z) — u{x)\ dx 



R2 



\u{xi + zi,X2 + Z2) - u{xi.,X2)\ da;idx2 
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which, by triangle inequahty, is less than or equal to 

\u{xi + zi,X2 + Z2) - u{xi,X2 + 2:2)1 dxida:2 

+ / \u{xi,X2 + Z2) — u{xi,X2)\ dxidx2 

< |zi| / \u{-,X2 + Z2)\bV(R) dX2 



+ \Z2\ / \u{xi,-)\BViR) d^l 

Jm 
< V2\z\\u\bv(r^), 
thanks to (IA.2I) and inequality \zi\ + \z2\ < \/2 |z|. D 
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